Simple Model of the Outflow from Starburst Galaxies: 
Application to Radio Observations 
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ABSTRACT 

We present an analytical model for the hydrodynamic outflow from the disk of a starburst 
galaxy. The model is used to calculate the cosmic ray propagation and the radio intensity 
distribution in the nuclear starbust region of NGC 253. We find that the cosmic ray energy 
production rate of the central 600 pc of NGC 253 is about 3 • 10 41 erg s _1 , that is about 15 
percent of the total mechanical supernova power. For this inner region we estimate a terminal 
outflow velocity of 900 km s _1 and a mass loss rate of M = 2 -j- AMp 



Subject headings: cosmic rays- acceleration- galaxies 

1. Introduction 

Observations of the soft X-ray and radio con- 
tinua and the H a line emission support the idea 
of galactic winds from starburst galaxies (see e.g. 
Strickland (2004) for a review). In these galax- 
ies which exibit a very high and spatially localised 
stellar formation rate, powerful stellar winds and 
supernova explosions produce a large amount of 
very hot gas. This gas expands out of the star- 
burst region and forms a supersonic galactic wind 
flow. 

Also in more quiescent galaxies such as the 
Milky Way galactic wind flows can arise. They are 
primarily driven by cosmic rays (Ipavich (1975); 
Breitschwerdt et al. (1987), (1991), (1993); Zi- 
rakashvili et al. (1996)). In such galaxies the gas 
heated by a spatially extended supernova activity 
loses its energy largely by radiative loses, whereas 
the nonthermal cosmic ray component, produced 
in the same context, can only cool adiabatically. 



VIq per year. 



In addition, the cosmic rays produced in the galac- 
tic disk can not freely escape from the Galaxy but 
rather amplify Alfven waves (e.g. Wentzel (1974)) 
that lead to an efficient dynamical coupling of the 
thermal gas to the energetic particle component 
(Ptuskin et al. (1997)). Therefore the cosmic ray 
pressure gradient is the dominant driver of such 
an extended galactic wind flow. Due to the differ- 
ent mechanisms these two types of outflow have 
somewhat different characteristics. 

Chevalier & Clegg ((1985)) introduced an ele- 
gant analytical hydrodynamic model for the galac- 
tic wind flow from a starburst galaxy. It assumes 
steady-state spherical symmetry. Despite its sim- 
plicity it seems qualitatively valid and is often used 
for estimates of the X-ray distribution and of the 
cosmic ray propagation. 

However, the assumed spherical shape of the 
starburst region is a disadvantage of this model. 
In many cases the starburst region has an ellip- 
soidal shape or a disk geometry. They are more 
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adequate for the problem considered. An analyt- 
ical approach to this case is not simple because 
one should solve at least a two-dimensional sys- 
tem of hydrodynamic equations. For this reason a 
number of numerical calculations for galactic wind 
flows have been made (e.g. Habe & Ikeuchi (1980); 
Tomisaka & Ikeuchi (1988); Suchkov ct al. (1994); 
Strickland & Stevens (2000)). They take also the 
temporal evolution of the flow into account. 

These studies evidently show that the power- 
ful energy release in the starburst region produces 
forward and backward shocks propagating in the 
galactic halo. They move mainly perpendicular 
to the galactic disk. After 3-10 million years 
(depending on the surrounding halo gas density) 
the shocks break through into the intergalactic 
medium. At smaller distances the galactic wind 
flow becomes steady state. It is bounded by al- 
most vertical walls of the shocked wind and halo 
gas (e.g. Suchkov et al. (1994)). 

The hydrodynamic equations can be signifi- 
cantly simplified for so-called self-similar station- 
ary flows. This means that different stream- 
lines have the same shape but different scales. 
The equations then reduce to ordinary differen- 
tial equations. This approach was first intro- 
duced by Bardeen & Berger ((1978)) for galactic 
winds. It has been successfully used in the the- 
ory of collimated rotating magnetohydrodynamic 
flows from accretion disks and of active galactic 
nuclei (e.g. Blandford & Payne (1984); Contopou- 
los & Lovelace (1994)). In spite of rather stringent 
boundary conditions at the base of the flow in the 
form of a spatial power-law gas pressure depen- 
dence, it appears to us as a useful tool. Indeed 
many nuclear starbursts do show power-law pres- 
sure profiles at distances of 100 to 1000 pc from the 
nucleus in the directions perpendicular and paral- 
lel to the major axis (e.g. Heckman et al. (1990)). 

In this paper we shall use such a self-similar 
approach for the description of the single-fluid hy- 
drodynamic outflow from a starburst galaxy. We 
have in particular found an explicit nontrivial an- 
alytic solution for one special (but not unreason- 
able) power-law dependence of the gas pressure at 
the base of the wind. 

We apply the solution found to a calculation 
of the radio emission from the central 600 pc of 
the starburst galaxy NGC253. The required ad- 
ditional transport equation for nonthermal rela- 



tivistic electrons couples the synchrotron emis- 
sion morphology, intensity, and spectrum with the 
overall outflow dynamics. This permits a determi- 
nation of the flow speed u and of the cosmic ray 
energy production rate Eqr as a measure of the to- 
tal mechanical energy production E$n ~ 10 Eqr 
in supernovae, through a comparison with the ra- 
dio synchrotron observations. We limit ourselves 
to this central region because a purely advective 
propagation model for the cosmic ray electrons is 
used. Such a model would not be adequate for 
the outer part of NGC253 since this outer region 
exhibits a steep radio spectrum (Hummel ct al. 
(1984), Carilli et al. (1992)) where diffusion dom- 
inates advection. Our approach may be considered 
as complementary to other methods related with 
X-ray and optical spectroscopic observations. 

The basic equations are given in the next sec- 
tion. A one-dimensional model for a thin starburst 
region similar to Chevalier & Clegg's model is de- 
scribed in section 3. Section 4 contains the ana- 
lytic solution which we found. Applications of this 
solution to cosmic ray propagation in the galaxy 
NGC 253 are given in sections 5 and 6. The out- 
flow parameters derived are described in section 
7. Section 8 contains a discussion of the results 
obtained together with our conclusions. 

2. Basic equations 

We shall assume that the starburst region lies in 
the galactic disk. The problem considered is mir- 
ror symmetric relative to the galactic midplanc. 
The galactic wind flow originates in the disk and 
extends into the galactic halo. In this region the 
steady-state axially symmetric gasdynamic equa- 
tions can be written in spherical coordinates R 
and 9 (the polar axis coinciding with the axis of 
galactic rotation): 
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Here ur and ug are the radial and latitudinal com- 
ponents of the gas velocity, respectively, 6 is the 
colatitudc, p and p are the gas mass density and 
pressure, respectively, and 7 is the adiabatic index 
of the gas. The first two equations are the Euler 
equations of motion in the radial and latitude di- 
rection respectively. Eq. (3) is the continuity equa- 
tion and Eq. (4) describes the evolution of the gas 
pressure. 

We neglect gravity, galactic rotation and dy- 
namical effects of the average magnetic field since 
they are insignificant for starburst outflows (cf. 
Chevalier & Clegg (1985)). Magnetic fields and 
rotation were taken into account in the framework 
of self-similar solutions for winds from accretion 
disks and active galactic nuclei. Gravitational po- 
tentials in a simplified form can also be included 
(cf. Bardeen & Berger (1978)). 

We shall seek the solutions of Eqs. (l)-(4) in 
the following form: 

u R {R,e) = u R {6)U{R/g(e)) (5) 
ug(R,6)=ue(6)U{R/g{6)) (6) 



p(R,e) = n(6)p (R/g(6)) 



(7) 



p(R,9)=n^9)P(R/g(6)) (8) 

Here the functions U(Rq), po(Ro), and P(Ro) de- 
scribe the spatial dependence of gas velocity, mass 
density and pressure at the base of the wind at 
6 = it/2, respectively, the function g{6) determines 
the shape of streamlines, whereas the functions 
ur(6), ug{6) and n{6) describe the dependence of 
the radial, latitudinal velocity components and gas 
density along a streamline, respectively. We put 
g(n/2) = 1 and n(7r/2) = 1. The equation for the 
streamlines is given by R = R g(6). The equation 
for the function g(6) may be written as 
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Substitution of Eqs. (5)-(8) into Eqs. (l)-(4) 
and taking into account Eq. (9) we obtain 
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Here ' denotes the derivative of the function P 
with respect to its argument i?o- Eq.(4) is reduced 
to an identity. From Eqs. (9) and (12) we find that 



ng ug sin 6 = const. 



(13) 



In order to separate the variables i? and 6 we have 
to assume that at the wind base P oc poU 2 oc Rq, 
with a constant S < 0. This means that the gas 
pressure at the base is a power-law function of 
radius and that the ratio of dynamical pressure 
to the thermal pressure (the square of the flow 
Mach number) at the base is constant. Without 
loss of generality we can put polI 2 /P = 1. Note 
that po(Ro) and U(Ro) need not individually fol- 
low a power-low in R , but rather only the product 
PoU 2 . Now Eqs. (10) and (11) are reduced to 
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Using Eqs. (9) and (13)-(15) one integral of the 
system can be found as 
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The Bernoulli constant for the each streamline 
differs from this integral by the factor U 2 (R/g(6j). 

As a result the problem is reduced to the so- 
lution of two ordinary differential equations (14) 
and (15) with the additional integral (16). If the 
solution is known, the function g can be found 
from Eq. (13). Boundary conditions should be 
determined at the wind base at 6 = ir/2. Gen- 
erally speaking, for arbitrary boundary values u R 
and ug, the solutions considered do not occupy 
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the whole space between 9 = and 9 = ir/2 (cf. 
Bardeen & Berger (1978)). The boundary values 
should be adjusted such that the solution occupies 
the whole space. 

As a natural boundary value we shall use ur = 
at 9 = 7r/2 for a starburst wind. It corresponds 
to a thin starburst region in the galactic disk. One 
can expect that the hot gas expands perpendicular 
to the disk in that case. 

3. Solution inside a thin starburst region 

If the starburst region has the form of thin disk, 
the one-dimensional hydrodynamic equations are 
a good approximation: 
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Here z is the height above disk midplane, q(z) and 
Q(z) are the sources of mass and energy, respec- 
tively. These equations are similar to the ones 
used by Chevalier & Clegg ((1985)) but for a plane 
parallel case. We shall assume that the problem is 
symmetric relative to the galactic midplane. 
Eqs. (17)-(19) can be integrated: 
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Here p(0) is the gas pressure in the midplane. 
The substitution of p and p from the first two 
equations into the third gives a quadratic equa- 
tion for the gas velocity u. Its solution is given 

by 
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The total energy and mass production rate per 
unit area of starburst disk are given by 



= 2 
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Here h is the halfwidth of the starburst disk. We 
can calculate the midplane gas pressure p(0) if the 
Mach number M = yj 'pu 2 Jyp at the boundary 
of the starburst disk z = h is given. Its value is 
determined by the wind solution in the halo region 
(cf. next Section). E.g. M = 1 in the Chevalier 
& Clegg model ((1985)). Using Eqs. (20)-(22) we 
found the gas pressure, the gas mass density and 
the gas velocity at the boundary of the starburst 
disk 
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These quantities formally coincide with the cor- 
responding quantities at the wind base in the 
thin disk approximation: p(h) = p(w/2), p{h) = 
p(n/2) and u(h) = —ug{-K/2). We should under- 
line that this solution is rather general and is not 
related to self-similarity of the halo outflow. It 
is valid for any distribution of energy and mass 
sources Qq{Rq) and qo{Ro) in the galactic disk. 
However, the value of M in general depends on 
Rq and is determined by a solution in the halo re- 
gion. For a starburst outflow M is close to unity if 
the adiabatic index is not close to unity (see also 
the next Section). 

4. Analytical solution in the galactic halo 

In the case 5 = —2j/(2j — 1) an analytical so- 
lution of Eqs. (14)-(16) can be found in the first 
quadrant < 9 < tt/2: 



UR = Moo cos6> 



(29) 
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u e = -(1 -7 _1 )u oc sin0 (30) 

n = sin 2 /^- 1 ) (31) 

g = sin-'"/^- 1 ) 6> (32) 

Here the value of Uqo fixes the asymptotic ve- 
locity at distances R that are large compared to 
the radius Rq for each streamline u ao U(R/ g(0)). 
The flow lines for this solution for 7 = 5/3 and 
5 = —10/7 and the surfaces of constant Mach 
number are shown in Fig.l in cylindrical coordi- 
nates r and z. The initial Mach number of the 
flow is M = ^2(7- l)/(2 7 - 1). For other val- 
ues of 5 the solutions of Eqs. (14)-(16) can be 
obtained numerically. The streamlines for these 
solutions look similar to lines shown in Fig.l. For 

5 closer to zero the streamlines are closer to the 
axis, i.e. the wind is collimated more strongly. We 
have not found the solutions for <5 < —27/ (27— 1) 
and for a natural boundary condition un ; = at 

6 = tt/2. Therefore the analytical solution with 
S = — 2 7 /(27— 1) has the maximally opened flow- 
lines in the framework of self-similar flows. 

As was first noted by Bardeen & Berger 
((1978)) the equations for self-similar flow do not 
contain any particular points corresponding to 
sonic points of the flow. The only particular point 
is the point where the sonic velocity equals the 
latitudinal component of the flow velocity (cf. Eq. 
(15)). This component is always smaller than the 
sonic velocity in the solutions found. The wind 
starts at the base with subsonic velocity (initial 
Mach number M w 0.756 for 7 = 5/3). It is accel- 
erated and passes through the sonic point at about 
9 « 70° (see Fig.l). This picture corresponds to 
a flow in a Laval nozzle. The cross-section of the 
Laval nozzle shrinks in the subsonic part of the 
flow in order to reach sonic velocities and expands 
after the sonic point in order to obtain an accel- 
eration of the flow beyond. As can be seen in 
Fig.l this is indeed the case. It is interesting that 
this flow geometry is organized by the gas itself 
through the perpendicular force balance. At the 
base of the wind the gas pushes the flowlines out 
of the axis. At larger heights the situation is re- 
versed. The wind originating closer to the axis 
expands faster. The gas pressure drops adiabati- 
cally and the gas pressure of the wind originating 
at larger distances pushes the flow lines to the 
axis. 




0.1 0.2 0.3 0.4 0.5 



Fig. 1. — Flow lines of the analytical solution 
{solid lines) for 7 = 5/3 and 5 = —10/7 and 
surfaces of constant Mach number (dashed lines). 
The units of the cylindrical coordinates r and z 
are arbitrary. The Mach number values are shown 
near the corresponding dashed lines. 

5. Propagation of cosmic ray electrons in 
the wind 

We shall neglect the diffusion of energetic elec- 
trons in the starburst outflow. This seems rea- 
sonable since the flow velocity is quite high. In 
addition, the cosmic ray scattering produced by 
self-excited Alfven waves (Ptuskin et al. (1997)) 
should be very effective in starbursts. This effect 
reduces the cosmic ray diffusion coefficient. The 
estimates of this effect in our Galaxy give for the 
diffusion coefficient the value 10 28 cm 2 s _1 for an 
energy of 10 GcV (cf. Ptuskin ct al. (1997)). It 
is inversely proportional to the cosmic ray pro- 
duction power per unit area of the galactic disk 
and is proportional to the magnetic field strength. 
Since starbursts have a high cosmic ray produc- 
tion power per unit area, typically a factor of 100 
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to 1000 or more larger than in our Galaxy, and 
since the magnetic field is a factor of 10 larger, 
the diffusion coefficient in starbursts is smaller by 
at least a factor of 10 to 100, and the diffusion 
is negligible in comparison with advection for 10 
GeV electrons which produce radio synchrotron 
emission ( for an application of these ideas to TeV 
gamma-ray emitting particles, see Aharonian et 
al. (2005)). The cosmic ray diffusion coefficient of 
much higher energy particles can in addition be- 
come quite small, if one takes into account that 
the star burst region is very turbulent, leading to 
additional strong magnetic field fluctuations which 
are produced by the turbulent velocity field. 

Assuming azimuthal symmetry and a steady 
state the isotropic part of the cosmic ray energy 
distribution function N(R, 9, E) obeys the follow- 
ing equation 
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Here the second term on the right hand side de- 
scribes the adiabatic energy losses, while b(E) de- 
scribes all other energy losses of particles. The 
boundary condition at the wind base 9 = ir/2 can 
be written as 



E 3 f) N f) 1 

(34) 

Here Nq(E) is the electron energy distribution at 
the wind base, bd(E) describes non-adiabatic en- 
ergy losses in the thin starburst region and qcn{E) 
is the differential source energy spectrum of cos- 
mic ray electrons normalized per unit area of the 
galactic disk. 

The synchrotron emissivity e u can be written 

as 



V3B 



dEN{E)f{v/v c 



(35) 



Here e and m are the charge and the mass of an 
electron respectively, B± is the strength of mag- 
netic field component perpendicular to the line 



of sight, v is the radio- frequency, c is the light 
velocity and / is the function that describes the 
electron synchrotron emission (cf. Landau & Lif- 
shitz (1975)). The characteristic frequency of syn- 
chrotron radiation v c is given by the formula 



3eB ± ( E 



Airmc \ mc 2 



(36) 



The synchrotron and inverse Compton energy 
losses can be written as 



m - ff) . (37) 

where ot = 87re 4 / (3m 2 c 4 ) is the Thompson cross- 
section, B is the magnetic field strength and U ra d 
is the radiation energy density. 

We shall neglect bremsstrahhmg and ionisation 
losses of electrons in the galactic wind flow. This 
is well justified, since the wind density is small. 
However, these losses may be non-negligible in the 
galactic disk which contains a large mass of gas. 
They are described by the second term on the left 
hand side of the boundary condition (34). The es- 
timates show that these losses can be significant 
in the very center of a starburst. This depends 
on the gas mass that is rather uncertain. In ad- 
dition it is not clear, whether radio electrons can 
well penetrate into gas clouds or whether they are 
quickly transported around clouds by the hot rar- 
efied wind. The cosmic ray penetration into the 
gas clouds in the starburst region is beyond the 
scope of this paper. For the sake of simplicity we 
shall neglect these energy losses in the disk also. 

The function / can be expressed in terms of 
the Mc-Donald function. A so-called (5-function 
approximation is often used for the function /. We 
shall use a more accurate method and approximate 
/ as (Pacholczyk (1970)): 



f(x) = 1.81a; 1/3 e" 



(38) 



where x = vjv c . 

The synchrotron intensity observed in a partic- 
ular direction is given by the integral of the emis- 
sivity along the line of sight: 
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The spectrum of sources of cosmic ray electrons 
is given by 



qcn(E) 



K-'QcRiRo) 



(40) 
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where Qcr is the total power that is transferee! 
by galactic sources to high energy cosmic ray pro- 
tons with energies between -E m i n and £ max per 
unit area of the galactic disk, and K is the pro- 
ton to electron ratio in the source. It is expected 
that such a spectrum is formed in supernova rem- 
nants by the diffusive shock acceleration mecha- 
nism (Krymsky (1977), Axford et al. (1977), Bell 
(1978), Blandford & Ostriker (1978)). 

We shall neglect the production of secondary 
electrons and positrons by the nucleonic cosmic 
ray component in the dense gas disk. Although 
this effect can be essential for starbursts (cf. 
Paglionc ct al. (1996)), its magnitude is rather 
uncertain because of the uncertain gas mass and 
the uncertain degree of penetration of the 10-^ 100 
GeV cosmic nucleons into the gas clouds. If signif- 
icant, this effect can be formally included in Eq. 
(40) by changing the proton to electron ratio K . 

The solution of equation (33) with velocity 
(5), (6) and energy losses b(R,6,E) = b (R,9)E 2 
can be written as 



N(R, 6, E) = N (E )^n 4 / 3 (6) (41) 

Here Eq = E m E / (E m — E) and the maximum 
energy E m of electrons is given by 
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Using the analytical solution (30) and (32) and 
the source spectrum (40) we obtain 
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where H{E) is the step function. The maximum 
energy E m is reduced to 
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6. Modeling of the radio spectra of NGC 
253 

NGC 253 is a nearby edge-on (inclination an- 
gle i = 78°; Pence (1981)) spiral starburst galaxy 
at distance d = 2.5 Mpc (Mauersberger et al. 
(1996)). It has extended radio (cf. Carilli et al. 
(1992)) and X-ray (cf. e.g. Strickland (2004)) ha- 
los that appear to be related with the powerful 
superwind flow originating in the galactic nucleus. 
Our model of the outflow can be used to model the 
spatial X-ray and radio distribution in the central 
region of this galaxy. 

We have performed the modeling of the radio 
distribution in the central 600 pc of NGC 253 us- 
ing the results obtained in the previous sections. 
The radio spectrum of this galaxy is flat at dis- 
tances Ro < 1 kpc from the nucleus (Carilli et al. 
(1992)) with a spectral index larger than -0.6. We 
expect that in this region cosmic rays propagate 
advectively. Concretely, we choose the function 
U(Rq) = 1. This means that the initial velocity is 
constant along the base of the wind and that the 
velocity depends on the colatitude 9 only. The an- 
alytical solution found in the previous sections can 
be written for this case as 



p(R, 9) = p g sin 
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(45) 



p(R,9) = _(2 7 -l)( 7 -l)pXooSin^# ( _^ 



2^ 
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Here the asymptotic velocity is the same for all 
streamlines, R g is the radius of the wind base and 
p g is the gas density at the wind base at R = R g . 
Dependence of velocity components on colatitude 
is given by Eqs. (29) and (30). 

We use a radiation field distribution at R < R g 
in the simplified form: 
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Here Lir and Lb are the infrared and starlight 
luminosities of the galaxy inside the region with 
wind base radius. It is reasonable to assume that 
radiation field is distributed similar to distribution 
of energy production in the starburst disk. 

For a calculation of the synchrotron emissiv- 
ity we need to define the magnetic field distribu- 
tion. Radio observations of NGC 253 (Beck et 
al. (1994)) and many other galaxies show that the 
radio emission is weakly polarized. This means 
that the random component of the magnetic field 
is larger then the regular component. In the case 
of a starburst this random component can be gen- 
erated by turbulence in the starburst disk and 
transported by the wind into the halo. Another 
possibility is to treat this random magnetic field 
as large amplitude Alfven waves propagating along 
the weak regular field and generated by the mo- 
tions of the footpoints of the field lines in the star- 
burst disk. This picture is similar to the case of 
the solar corona. 

We shall treat the random magnetic field as a 
fluid with an adiabatic index j m . It is known that 
Alfven waves have the adiabatic index j m = 3/2 
(Dewar (1970), McKee & Zweibel (1995)); an 
isotropic random magnetic field has the adiabatic 
index j m = 4/3 (Mestel (1965)). We used the fol- 
lowing expression for the magnetic field which has 
a similar R-dependence a,sp(R,9) (cf. Eq. (46)): 



B = B (sin 6) 7 



(48) 



Here B is the magnetic field strength at R = R g . 
The magnetic field component B± in Eqs. (35) 
and (36) was calculated as B± = y/2/'SB. 

The distribution of cosmic ray sources in the 
starburst disk in Eq. (43) was taken in a similar 
form 
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where Lqr is the cosmic ray energy production 
rate of the starburst. 



The continuum radio emission of galaxies con- 
sists of two main components: synchrotron emis- 
sion of cosmic ray electrons and a thermal emis- 
sion of the ionized gas which is produced mainly 
in HII regions, ionized by OB stars (see e.g. Con- 
don (1992) for a review). The thermal emission is 
closely related with thermal (so-called free-free) 
absorption which can be significant at low fre- 
quencies and can produce a flattening and even 
a turnover of the spectrum. It appears that this 
indeed the case for the spectrum of the nucleus of 
NGC 253 (Carilli (1996)). Other possible reasons 
for the flattening are reacceleration of low-energy 
electrons (Pohl & Schlickeiser (1990)), and ionisa- 
tion and brcmsstrahlung energy losses of electrons 
in the starburst disk (Hummel, (1991)). We give 
preference to free-free absorption here, since the 
thermal component of the radio emission is un- 
doubtedly present in NGC 253. 

To take free-free absorption into account we 
shall assume that the thin galactic disk seen face- 
on has the following distribution of the absorption 
depth 



2-1 
" 2 7 -l 



1 GH2 



(50) 



Here t is the absorption depth of the disk at 1 
GHz and at Ro = R g - It was assumed that this 
distribution is similar to the distribution of sources 
of energy in the galactic disk. This is reasonable, 
since the absorbing gas is ionized by OB stars. The 
number of these stars is proportional to the input 
of energy. The inclined disk of NGC 253 absorbs 
the synchrotron emission of the northern part of 
the halo. 

The thermal emission brightness of the inclined 
disk can be expressed in terms of the absorption 
depth r and the inclination angle i: 



Jt\i 



Tv 1 



1 



T j \ cos i I 



)• 



(51) 



Here T is the absorbing gas temperature in en- 
ergetic units. Generally speaking, the two last 
equations are valid if the absorbing gas has a large 
covering factor. This is possibly not the case in the 
disk of NGC 253. However it is possible in the cen- 
tral region, where the main part of the absorption 
occurs. 

We should take into account the synchrotron 
emission of supernova remnants (SNRs) in the 
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starburst disk. The SNRs in the disk have a 
small volume filling factor. However, the mag- 
netic field strength and the number density of rela- 
tivistic electrons in SNRs are larger in comparison 
with the corresponding quantities in the interstel- 
lar medium of the starburst disk and in the wind 
halo. We used the following expression for the 
brightness of the radio emission produced by the 
SNRs from the inclined disk: 

T _ 7-1 Jod 2 (R \~^ ( v N-0-5 

J SNR- 27 _ l7r ^ r {jfj {TgrJ 

x (l -e- r/|cosi| ) . (52) 

This emission from SNRs is also influenced by free- 
free absorption. The total radio flux J from SNRs 
is about one tenth of the total radio flux for normal 
galaxies (Condon (1992)), but it can be larger in 
galaxies with a fast escape of electrons. A more 
quantitative estimate was obtained by Ulvestad 
((1982)) using the relation between radio surface 
brightness and diameter (the "£ — £)" relation, see 
also Bcrezhko & Volk (2004)): 

v / _BsN \ ~ 1/17 / niSM \ - 2 / 17 , . 
H^ergJ llcm-3j [M) 

where z/sn is the supernova rate. This flux weakly 
depends on the interstellar medium number den- 
sity nisM and the supernova explosion energy i?sN- 
We recalculated the numerical factor in Eq. (53) 
using the X— D relation for supernovae in the star- 
burst galaxy M82 (Huang et al. (1994)). It is a 
factor of two larger than the factor obtained by 
Ulvestad ((1982)). That means that supernovae 
in starburst galaxies are slightly brighter than in 
normal galaxies. This number is in a good agree- 
ment with radio observations of galaxies (cf. Con- 
don (1992) for a review) and theoretical estimates 
(Lisenfeld & Volk (2000)). 

We have performed the modeling of the star- 
burst wind with a radius R g = 300 pc and have 
calculated synchrotron and thermal emission using 
the formulae given in this and previous sections. 
The choice of this radius was done in accordance 
with the size of the region where the radio emission 
of NGC253 is concentrated (Ulvestad (2000)). 



For the adiabatic index of the random mag- 
netic field we use the adiabatic index of Alfven 
waves j m = 3/2. Since cosmic ray and Alfven 
wave pressures are not negligible (see below) we 
also used a gas adiabatic index 7 = 3/2. The tem- 
perature of the absorbing gas was fixed at 8 • 10 3 
K. The infrared luminosity of the nucleus of NGC 
253 was taken as Lir = 1.24 • 10 10 L Q (Radovich 
ct al. (2001)). For the starlight luminosity, we 
used the luminosity of the bulge of NGC 253, 
L B = 0.39 • 10 10 L© (Pence (1980)). Also a source 
proton-to-electron ratio K = 70 was chosen. Fi- 
nally, a maximum energy B miK = 10 6 GeV of the 
cosmic ray protons was assumed. The results are 
shown in Fig. 2 and Fig. 3. 

The contour plot of the 1.4 GHz brightness dis- 
tribution modeled is shown in Fig. 2. The radio 
distribution was convolved with a Gaussian beam 
with FWHM 5.1" x 3.8" (62 x 46 pc at 2.5 Mpc 
distance), and an inclination angle of i = 78° was 
used. The distribution is slightly asymmetric rela- 
tive to the horizontal axis due to free-free absorp- 
tion. The comparison with the observed distribu- 
tion (Mohan et al. (2005)) permits us to estimate 
the wind velocity. We found that the wind velocity 
■Moo should be larger than 500 km s _1 . For smaller 
velocities the brightness drops too fast with height 
above midplane because of strong synchrotron and 
inverse Compton losses. We have found a satisfac- 
tory fit for L CR = 2.94 • 10 41 erg s _1 , = 900 
km s _1 , r = 0.05| cosi|, J = 1.0 Jy and B = 50 
^G. The magnetic field pressure was taken close 
to the cosmic ray pressure. The value of the ra- 
dio flux Jo from SNRs corresponds to the super- 
nova rate z^sn = 0.06 yr^ 1 which is intermediate 
between the value z^sn = 0.03 yr _1 of Matila & 
Meikle ((2001)) and z/ S n = 0.1 yr" 1 (Ulvestad & 
Antonucci (1997)). Finally, a mean number den- 
sity nisM = 100 cm~ 3 of the interstellar medium 
in NGC253 (Paglione et al. (1995)), and a mean 
energy of supernova explosions i?sN = 10 51 erg 
were assumed in Eq. (53). 

The uncertainty of the velocity derived in this 
way is larger from the side of higher velocities. 
For example, for u^, = 500 km s _1 the contour 
lines near the upper and lower parts of the plot are 
shifted by about one contour level in the midplane 
direction. A similar shift in the opposite direction 
corresponds to = 2000 km s _1 . The derived 
value — 900 km s _1 produces a good simi- 
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larity between observed and modeled radio distri- 
bution in the southern part of the outflow. We 
think that the radio distribution in the northern 
part of the outflow should be explained by using a 
more complex model with slightly larger and vari- 
able wind velocities along the base (for example, 
another function U(Ro) might be used). This can 
produce spur-like features in the observed radio 
distribution in the northern part of the outflow. 
Another possible explanation is limb-brightening 
due to the superwind's lateral interaction with 
halo gas which is not included in our model. 

The spectrum of the NGC 253 nuclear region is 
shown in Fig. 3 (solid line). The data points with 
the same resolution 33" x 21" (400 x 255 pc at 2.5 
Mpc distance) are taken from Carilli ((1996)). 

The total radio flux of SNRs amounts to about 
one third of the total radio flux at 1.4 GHz. The 
thermal emission amounts to one third of the total 
intensity at 10 GHz. This is in good agreement 
with results for the separation of the thermal and 
nonthermal emission for the whole galaxy NGC 
253 (Niklas et al. (1997)). The fit can be better if 
a smaller temperature of the emitting gas is used. 
It should be noted that Israel & Mahoney ((1990)) 
found that the flattening of the radio spectra of 
galaxies at small frequencies can be explained by 
free-free absorption in a hypothetical cold ionized 
gas. 

The observed synchrotron intensity mainly con- 
strains the ratio of the cosmic ray power Lqr to 
the wind velocity. The wind velocity and cosmic 
ray power influence the spectral shape of the radio 
emission only weakly. 

The cosmic ray luminosity found corresponds 
to the cosmic ray pressure at the wind base (at 
R = R g ): 



PGR 



7 



L 



CR 



8(2 7 -l)7r Moo i?2 



(54) 



This gives pgr = 1-08- 10~ 10 erg cm~ 3 at the wind 
base at Ro — R g . The total nonthermal (cosmic 
ray + Alfven waves) luminosity is about 5.6 • 10 41 
erg s _1 . 

7. Estimates of the outflow parameters 

We have modificated the Chevalier & Clegg 
model ((1985)) for a thin starburst region (cf. 
Sect. 3). In the particular case of the analytical 
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Fig. 2. — Radio brightness distribution of NGC 
253 at 1.4 GHz calculated in our model. The mod- 
eled radio intensity was convolved with a Gaussian 
beam of FWHM 5.1" x 3.8" (62 x 46 pc at 2.5 Mpc 
distance) . The solid contour levels are a geometric 
progression in 1.7, ending at 0.446 Jy/beam in the 
center. The projected distance in pc is shown on 
both axes. The observed brightness distribution 
at 1.4 GHz (Mohan et al. (2005)) is shown by the 
gray scale contours. 



solution (29), (30), (45), (46), the Eqs. (26) - (28) 
can be rewritten as 
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7 



E 



■^juih) = 5.63- 10 8 cm s" 1 



1/2 



I 10 43 erg s" 1 I \ M Q yr 
p{h) = 1.88-10~ 27 g cm 
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(56) 



p(h) = 2.98-10" 10 erg cm" 37 1 



E 



1/2 



7 2 \ 10 43 erg s" 1 



M 
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M Q yr- 



R„ 



Ro\ 
1 kpc ) { R c 



(57) 

Here E and M are the total energy and mass pro- 
duction rates of the starburst disk at R < R g . 



10 



2 


J(l>), Jy/beam 




1 




I 


0.5 










v , GHz . 


0.1 1 


10 



Fig. 3. — Spectrum of the nuclear region of NGC 
253 (solid curve), calculated in our model. The 
modeled radio intensity was convolved with a 
Gaussian beam of FWHM 33" x 21" (400 x 255 pc 
at 2.5 Mpc distance). The data points are taken 
from Carilli ((1996)). 

The energy production rate is determined 
mainly by supernovae. The upper limit for the 
supernova rate in the nucleus of NGC 253 is 0.3 
yr _1 (Ulvestad & Antonucci (1997)). It seems 
that the actual number is not larger than 0.1 yr" 1 
because of the absence of a new supernova during 
the last years. For a supernova rate 0.06 yr _1 in 
NGC 253 the energy production rate E = 1.9- 10 42 
erg s _1 under the assumption of negligible radia- 
tive losses of the gas. The nonthermal luminosity 
found in the previous section is about one quarter 
of this number. 

Since part of the supernova energy can be lost 
radiatively, the energy input may be smaller. The 
X-ray observations give an the energy input E w 
10 42 erg s" 1 (Strickland et al. (2002)) and an av- 
erage wind velocity u between 400 and 500 km 
s _1 that is between the velocity 300 km s _1 at the 
wind base and the asymptotic velocity 900 km s _1 
found in the previous section. 

For this case we found the mass-loss rate M = 
3.9M Q yr -1 , and the total pressure p = 4.6 • 10~ 10 
erg cm~ 3 at R = R g and a factor of 15 larger at 
a distance of 50 pc from the nucleus. These num- 
bers should be compared with the cosmic ray and 
magnetic pressures derived in the previous section 
p C R = 1.1-10" 10 erg cm" 3 and B 2 /87r = LO-IO" 10 
erg cm -3 respectively. Therefore the nonthermal 
pressure is almost one half of the total pressure. 
The rest can be attributed to the thermal gas pres- 



sure, that is p g = 2.5 • 10~ 10 erg cm~ 3 . The 
corresponding mass density is p g = 3.8 • 10~ 25 g 
cm~ 3 . The corresponding temperature of a one- 
phase fully ionized gas is about 5 • 10 6 K. 

This is not surprising. If 15 procent of the su- 
pernova energy is converted into cosmic ray en- 
ergy and the rest is partially radiated, then the 
gas pressure may be comparable with the cosmic 
ray pressure. In this case it is reasonable to use 
a smaller adiabatic index 7 < 5/3 of the gas to 
include cosmic rays and magnetic fields. 

In reality, a one-phase thermal gas with a tem- 
perature of several million degrees is excluded by 
the X-ray observations, since it would produce 
too large an amount of X-rays (Strickland et al. 
(2002)). The main fraction of the thermal en- 
ergy may rather be contained in an even hotter 
gas component, and the main part of the mass 
may then be contained in gas with a temperature 
below one million degrees. 

We can even consider an extreme case by as- 
suming that the gas radiates all its energy. The 
superwind flow is then driven by cosmic rays and 
Alfven waves (or turbulent magnetic fields) alone. 
Using a nonthermal luminosity E = 5.6 • 10 41 erg 
s _1 , as found in the last section, we then obtain a 
mass-loss rate of M = 2.2 M yr _1 . 

8. Conclusion 

We have presented a simple analytical model 
for a steady state axially symmetric gas-dynamic 
outflow from a starburst galaxy. In the case of a 
power-law radial distribution of the gas pressure at 
the wind base self-similar solutions can be found. 
For the special choice of the power-law index S = 
—27/(27— 1) « 1.5 we found an explicit analytical 
solution. This value of the index seems reasonable, 
since many starbursts show a strong concentration 
of the energy release towards the center. 

The solution found has less steep pressure and 
density profiles in comparison with the Chevalier 
& Clegg model (Chevalier & Clegg (1985)). It is 
well known that the Chevalier & Clegg solution 
is in agreement with the X-ray brightness of M82 
at distances larger than 1-2 kpc from the nucleus, 
whereas at smaller distances the observed bright- 
ness relief is lower (e.g. Fabbiano (1988)). This 
may be understood in the framework of our model. 
At distances smaller than the starburst region gas 
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collimation is possible. This results in weaker de- 
pendencies on radius. At larger distances the wind 
geometry is roughly spherically symmetric, per- 
pendicular forces play no role and the Chevalier & 
Clegg model is more adequate. 

In reality the power-law pressure dependence 
can not be valid at small radii either, i.e. there 
exists some minimum radius. The wind originat- 
ing inside this radius is not described by the self- 
similar solution. However, for a small minimum 
radius it blows in a narrow region near the axis. 
Therefore this restriction appears not to be essen- 
tial. 

The solution found was applied to the model- 
ing of the radio distribution of the central part 
(R < 300 pc) of the starburst galaxy NGC 253. 
This galaxy contains a very bright radio nucleus 
with a flat spectrum (with spectral index rnaging 
from —0.43 to —0.17 at different frequencies) and 
a less bright extended radio disk with a steeper 
spectrum of spectral index -0.7 (Hummel et al. 
(1984), Carilli et al. (1992)). 

The propagation of the radio electrons can be 
purely advective in the central part of NGC 253 
since the gas velocity is rather high (about 1000 
km s" 1 ). It appears that this velocity is smaller 
in the outer part of the galactic disk (R > 1 kpc) 
where the spectrum is steeper. It is possible that 
the outer part contains a cosmic ray driven galac- 
tic wind which has small velocity at the wind base 
(cf. Breitschwerdt et. al. (1991); Zirakashvili et 
al. (1996)) and where the diffusion of particles 
results in the spectral steepening. 

Our model is applicable for the central super- 
wind region (R < 300 pc) of NGC253 which con- 
tains the central starburst. NGC253 has also 
a comparable energy release at larger distances. 
However this is distributed over an area which is 
at least a factor 100-1000 larger than the central 
starburst. Our model does not describe that re- 
gion which contains a quasi-static galactic halo or 
a slow cosmic-ray driven galactic wind. It also 
does not describe the interaction of the gas in this 
region with the superwind flow. 

We have obtained a satisfactory fit for the radio 
brightness distribution and for the radio spectrum 
of the nucleus of NGC 253. This permits us to 
estimate the parameters of the outflow: the wind 
velocity and mass loss rate of the starburst. 



From a physics point of view these quantities 
should be determined from a theoretical model 
which describes how thermalised supernova ejecta 
and the interstellar gas heated by supernova 
shocks radiate some part of the energy and form 
an outflow from the starburst disk. In the absence 
of such a theory we limit ourself to the determi- 
nation of these quantities from observations. 

We found an asymptotic velocity of 900 km s _1 
from this inner region. For an energy input of 
10 42 erg s _1 in the starburst nucleus of NGC253 
we obtain a mass loss rate of 4M Q yr _1 . This 
number drops down to 2M yr _1 for the extreme 
case when the gas radiates almost all its energy. 
Even this case does not contradict the X-ray obser- 
vations, since it appears that the X-rays are pro- 
duced by a hot gas that is formed from the interac- 
tion between the superwind and the circumstellar 
medium (halo gas, clouds etc.) (cf. Strickland et 
al. (2002)). 

The cosmic ray production power is estimated 
as 15 procent of the supernova power. We found 
that cosmic ray and magnetic pressures are com- 
parable with the gas pressure. Therefore we can 
conclude that the starburst outflow in NGC 253 
is collectively driven by hot gas, cosmic rays and 
turbulent magnetic fields. 

We thank the anonymous referee for a number 
of valuable suggestions regarding the interpreta- 
tion of the theoretical results. 
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